The Quark Propagator in the NJL Model in a self-consistent l/N c Expansion 
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The quark propagator is calculated in the Nambu-Jona-Lasinio (NJL) model in a self-consistent 
l/iV c -expansion at next-to-leading order. The calculations are carried out iteratively in Euclidean 
space. The chiral quark condensate and its dependence on temperature and chemical potential 
is calculated directly and compared with the mean-field results. In the chiral limit, we find a 
second-order phase transition at finite temperature and zero chemical potential, in agreement with 
universality arguments. At zero temperature and finite chemical potential, the phase transition 
is first order. In comparison with the mean-field results, the critical temperature and chemical 
potential are slightly reduced. We determine spectral functions from the Euclidean propagators by 
employing the Maximum-Entropy- Method (MEM). Thereby quark and meson masses are estimated 
and decay channels identified. For testing this method, we also apply it to evaluate perturbative 
spectral functions, which can be calculated directly in Minkowski space. In most cases we find that 
MEM is able to reproduce the rough features of the spectral functions, but not the details. 



I. INTRODUCTION 

Describing the low-energy sector of QCD is a com- 
plicated task due to its strong coupling. Especially in- 
vestigating the region of the chiral and deconfinement 
phase transition requires non-perturbative methods. To 
describe these phenomena, one therefore often employs 
effective models, which give a simplified description of 
the theory and are valid in a limited energy range. In 
this article, we focus on the Nambu-Jona-Lasinio (NJL) 
model [l| , where the quark-gluon interactions of QCD are 
substituted by effective 4-quark interactions. This model 
does not exhibit confinement, which is its major lack, 
but it incorporates chiral symmetry. Spontaneous and 
explicit breaking of this symmetry as well as its restora- 
tion at high temperatures or densities can be realized. 

Despite the simplifications in the interaction, the NJL 
model cannot be solved exactly, but further approxima- 
tions are necessary!]] In most publications, the model has 
been treated in mean-field (Hartree) and random-phase 
approximation to describe the chiral dynamics of quarks 
and the meson spectrum, both, in vacuum and in hot and 
dense matter [3H5| . A shortcoming of these approxima- 
tions is that the effects of mesons on the quark propagator 
are not included. In the quark sector this leads to the 
wrong universal behavior at the phase transition and pro- 
duces delta peaked quark spectral functions. As a con- 
sequence, for instance, they are not suitable to be used 
in the Kubo formula for the calculation of shear viscosi- 
ties Q. In the mesonic sector, hadronic decay channels, 
like p — > 7T7T, are not included, which are the physical 
decay channels in the confined phase. 



1 Strictly speaking, since the NJL model is non-renormalizable, a 
unique exact solution does not even exist, but the results depend 
on the regularization scheme. An interesting alternative to the 
continuum methods discussed in this article is to solve the NJL 
model on the lattice Q. In this case the non-renormalizability 
has the consequence that there is no continuum limit, so that 
the results depend on the choice of the lattice. 



These processes can be included systematically within 
a l/iV c -expansion, i.e., within an expansion in the inverse 
number of colors, beyond the leading order, see, e.g., 
Refs. [7-16|. Here one can basically follow two different 
approaches: In the "strict 1/N C expansion scheme" one 
first solves the gap equation in leading order, i.e., Hartree 
approximation, and then adds the 1 /N c corrections to the 
quantities of interest perturbatively, without modifying 
and solving the gap equation again. This method yields 
good results, e.g., for the rho- meson dominated electro- 
magnetic pion form factor in the time-like region [lij as 
well as for the low-temperature behavior of the quark 
condensate [ill EH an d the pressure Above the 

critical temperature the model has been studied in [I?} 
where the influence of soft modes on the quark spectral 
function was investigated. On the other hand, the pertur- 
bative approach breaks down in the vicinity of the phase 
transition, where a method is needed which incorporates 
the 1/N C corrections self-consistently in the gap equation. 
First attempts in this direction have been performed in 
Refs. [il dHHH , but in a simplified approach, where non- 
local contributions to the quark self-energy have been ne- 
glected. This scheme is thermodynamically inconsistent 
and it was found that the chiral phase transition at finite 
temperature is first order for two quark flavors (13 . [l5j , 
in contradiction to universality arguments [18j. 

In the present paper, we derive a self-consistent solu- 
tion of the gap equation using the so-called 4>-derivable 
theory in next-to-leading order without further approxi- 
mations. This approach is thermodynamically consistent 
and allows a meaningful investigation of the phase tran- 
sition. Due to the nonlocal self-energy contributions, the 
structure of the equations is rather involved and we use 
the imaginary time (Matsubara) formalism to keep them 
on a tractable level. This allows us to study static prop- 
erties, like the quark condensate at zero and finite tem- 
perature. On the other hand, the analytic continuation 
of dynamical quantities to real times is problematic. In 
order to gain information about quark and meson spec- 
tral functions in Minkowski space, we therefore employ 



the Maximum- Entropy Method (MEM). 

The remainder of this article is organized as follows. 
In Sec. Hi] we summarize the basics of the model and the 
leading-order formalism. In Sec. IIIII we introduce the <fr- 
derivable theory and formally derive the gap equation for 
the quark propagator in next-to-leading order (NLO) in 
1/N C . The numerical results for the self-consistent solu- 
tions in the Matsubara formalism are shown in Sec. IIVI 
There, we also discuss the behavior of the quark conden- 
sate as a function of temperature and chemical potential. 
In Sec. [V] we review the basic ideas of MEM and apply 
them to study quark and meson spectral functions. We 
conclude with a summary in Sec. IVI1 

II. MODEL AND LEADING ORDER 
FORMALISM 

We use a two-flavor NJL model with a scalar and pseu- 
doscalar interaction, given by the Lagrangian 

C = q{i$ - m )q + G[(qqf + {qi lb fqf] (1) 

with a dimensionful coupling constant G and the Pauli 
matrices f in isospin space, mo is a small bare quark 
mass, which explicitly breaks chiral symmetry. For calcu- 
lations in the chiral limit it is set to zero. In this limit the 
Lagrangian is invariant under SU(2)l x SU(2)r transfor- 
mations. For nonvanishing but small values of mo, this 
is still an approximate symmetry of the model. 

The 1/N C counting scheme is introduced in the NJL 
model by assuming that the quark fields q have N c color 
degrees of freedom. Consequently, a closed quark loop 
yields a factor N c . Furthermore, it is assumed that the 
coupling constant G scales like 1/N C . In this article these 
rules are used to organize the diagrams in a systematic 
and symmetry conserving way. In all explicit calcula- 
tions, however, we take the physical number of colors, 
iV c = 3. 

In any approximation the full quark propagator S(k) 
is given by 

S-\k) = Sv\k)-nk) (2) 

with the inverse bare propagator Sq (k) — # — mo and a 
self-energy For large enough couplings, chiral sym- 

metry is spontaneously broken and the quarks acquire 
a dynamical mass much larger than the bare mass. In 
leading order in 1/N C this is described by the Hartree 
Dyson-equation, Fig. [T] In the self-energy insertion on 
the right-hand side (r.h.s.), the factor N c of the quark 
loop is compensated by a factor 1/N C from the vertex. 
Hence, if we assume that the bare propagator is of the 
order N® , we find that the dressed propagator in Hartree 
approximation is strictly of the order N® as well. Note 
that in this paper we draw the local four-point interac- 
tion as a wavy line, indicating the direction of the inter- 
action. In the present example this prevents confusion of 
the Hartree term with the Fock term, which is suppressed 
by one order of l/N c . 
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FIG. 1. Dyson-equation for the Hartree quark propagator 
(bold lines). Thin lines represent the bare propagator, the 
wavy line the bare interaction. 

In vacuum, the Hartree self-energy is given by 

E ff = 2zG J -^Tr(S(k)) . (3) 

It is local and purely scalar and therefore corresponds to 
a constant shift in the quark mass, 

n%H = m + Sff . (4) 

The dressed or "constituent quark" mass m# is the scalar 
part of the inverse Hartree propagator S^ 1 (k) = $ — 
rriH- In this way it enters the r.h.s. of Eq. ([3]), thereby 
giving rise to a self-consistency problem. 

To describe the system at nonvanishing temperature T 
and chemical potential /i, we apply the Matsubara for- 
malism. The quark propagator is then defined at discrete 
imaginary energies iw n + fj,, with fermionic Matsubara 
frequencies uj n = (2n + 1)ttT. Accordingly, the energy 
integration in fermionic loop integrals is replaced by a 
sum, 

In the Hartree approximation, where the self-energy is 
just a constant, the analytic continuation of the propa- 
gator to real energies is, of course, trivial. However, this 
will no longer be the case at NLO. 

Mesons are described by a Bethe-Salpeter-equation for 
the quark- antiquark T-matrix, 

if(q) = iK + iK{-ih(q))iT(q) . (6) 

The leading order corresponds to the random-phase ap- 
proximation (RPA) , depicted in Fig. [2J In this case 

%k = 2iG^(r M ®r M ) (7) 

M 

is the bare scattering kernel and 
Jm(q) =T M Il(q)T M 

= i I ^Tv[T M S(k + q)T M S(k)} (8) 




FIG. 2. Bethe-Salpeter-equation for quark-antiquark scatter- 
ing. Double lines denote RPA meson propagators. 



are the quark-antiquark polarization functions in the 
scalar and pseudoscalar channels, 



T s = 1, T P:a = ij 5 r a . 



(9) 



Here we have used that these channels do not mix. 
Eq. ^ is then easily solved with the ansatz 



T M = -D m (Tm®T m ) 



which yields 



D M (q) = 



2G 



1 - 2GJ M (q) 



(10) 



(11) 



Following Refs. [14|,|l5j], we will call the functions -Dm(<7) 
"meson propagators" although they are not properly nor- 
malized. In particular, the meson masses are given by the 
pole positions of the propagators, 



= 0. 



(12) 



When we expand the r.h.s. of Fig. [5] into a geometric 
series, we see that the term with n quark loops contains 
ri + 1 four-point vertices. Thus, applying the 1/N C count- 
ing rules, we find that the RPA meson propagators are 
strictly of the order 1/N C . 

In the Matsubara formalism, the meson propagators 
and polarization functions are defined at discrete imag- 
inary energies iui m , with bosonic Matsubara frequencies 
Lu rn = 2nnrT. The determination of meson masses ac- 
cording to Eq. (fT2"j) then requires the analytic continu- 
ation of the propagator to real energies. Again, this is 
easily done in the Hartree + RPA scheme, but will be 
non-trivial at NLO. 



III. ^-DERIVABLE THEORY 

As motivated in the introduction we are aiming at a 
self-consistent extension of the approximation scheme be- 
yond the leading order in 1/N C . To this end we apply 
the l/iV c -expansion on the level of the thermodynamic 
potential using the so-called ^-derivable theory 0, [2(| ■ 
This scheme preserves all important symmetries. 

The full thermodynamic potential is given by 

0[5] = iTrlnS- 1 +iTr{Y,S) + $[S}, (13) 



where S and E = Sq — S^ 1 are the full quark propagator 
and the full self-energy, respectively, and Tr denotes a 
functional trace over all space-time and internal degrees 
of freedom. The functional $[S I ] summarizes all closed 
two- particle irreducible (2PI) diagrams plj . 

The stationarity condition -ML. = implies that 



E(z) 



5{iS(x)) 



(14) 



i.e., the self-energy can be obtained as a functional 
derivative of $. Diagrammatically, this corresponds to 




FIG. 3. Leading order contribution to <E>. 



cutting a single quark line of $ at all possible places. 
In turn, $ depends on E via the full quark propagator. 
Eq. (fT4|) therefore constitutes a self-consistency problem. 

Similarly, the symmetry conserving scattering kernel 
for the mesonic BSE can be obtained as 



K(x,y) = - 



5 2 <S> 



6(iS(x))6(iS(y)) ' 



(15) 



corresponding to cutting the ^-functional twice. How- 
ever, unlike for the self-energy, this is not a self- 
consistency problem because $ does not depend on K. 

Self-consistent approximation schemes can now be in- 
troduced by performing truncations of the functional $. 
In the present article, we expand $ in powers of 1/N C 
to next-to-leading order. The leading-order contribution 
corresponds to the "glasses" diagram, shown in Fig. [3] 
In vacuum it is given by 



M ^ J 



rf 4 fc 



Tr(r M 5(fc)) 



(16) 



which can be generalized in the Matsubara formalism by 
the replacement (JSJ) - As it contains two quark loops and 
one four-point vertex, it is of the order 0(N C ). Cutting 
one or two quark lines, we reproduce our earlier result 
that the Hartree self-energy and the RPA scattering ker- 
nel, respectively, are the corresponding leading-order ex- 
pressions in this expansion. Formally, this can also be 
obtained from Eqs. (fT4|) and (fl5|) . where one has to take 
into account that the transformation to momentum space 
brings in extra factors of (27r) 4 . For instance, for the self- 
energy one gets 



E(jfe) = i(27r) 4 



S<S> 

SS(k) ' 



(17) 



which, when applied to Eq. (|T^|) . indeed yields Eq. 

The fact that the Hartree self-energy and the RPA 
scattering kernel can consistently be derived from the 
same ^-functional guarantees that this approximation 
scheme is symmetry conserving. In particular chiral 
Ward identities and low-energy theorems are fulfilled and 
the RPA pion is massless in the chiral limit, as required 
by Goldstone's theorem. 

The NLO correction to $ is given by the ring sum, 
depicted in Fig. |4] Taking into account the appropriate 
symmetry factors, these diagrams can be combined to a 
logarithm. One obtains 



M J 



d 4 q 



ln(l - 2GJ M (q)) , (18) 
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FIG. 4. NLO contribution to $. 



which depends on S through the polarization function 
J Mi Eq. ©. Applying again Eq. (IT71) . we find the fol- 
lowing NLO correction to the self-energy: 

E«(fe) = -^D M (q)T M S(k-q)r M . (19) 

It describes the dressing of the quark propagator by an 
RPA meson and corresponds to the insertion in the last 
diagram in Fig. [5] Recalling that the RPA meson propa- 
gators are of the order l/N c , this self-energy term yields 
a correction of the order 1/N C to the quark propagator. 
However, when the diagrams are iterated in the gap equa- 
tion, as shown in the figure, higher orders are generated. 
Therefore both, the self-consistent quark propagator and 
the individual self-energy contributions, are no longer of 
strict orders in 1/N C . Obviously, the same is true for 
and which contain higher orders in l/N c as well, 
when the self-consistent quark propagator is used. In the 
present scheme, the 1/N C counting is thus introduced on 
the level of skeleton diagrams for the ^-functional, i.e., 
before dressing the propagators. 

As we have discussed, at leading order RPA mesons to- 
gether with the Hartree gap equation are consistent with 
chiral Ward identities so that, in the chiral limit, pions 
emerge as massless Goldstone bosons in that scheme. In 
the same way one can construct mesons which are consis- 
tent with the NLO gap equation. To that end one has to 
calculate the NLO corrections to the scattering kernel by 
applying Eq. (TT51) to corresponding to cutting the 

diagrams in Fig. [4] twice, and iterate them together with 
the leading-order kernel in the BSE. Again, this scheme 
preserves chiral symmetry and the resulting pions are 
massless in the chiral limit. 

In this context it should be noted that the mesons 
which enter the NLO self-energy diagram are not the 
NLO-corrected mesons, but RPA mesons. To be precise, 
they are obtained from the BSE with the leading-order 
scattering kernel, Fig. [2j but the polarization functions 
J Mi Eq. ©, involve the self-consistent NLO quark prop- 
agators. As a consequence, these "intermediate mesons" 
are not restricted by chiral Ward identities and the pi- 

- + -9 + c\ 

FIG. 5. NLO gap equation. The thin and bold lines represent 
bare and dressed quark propagators, respectively. The double 
line symbolizes RPA- like mesons as defined in Fig. [2] but in- 
volving the self-consistent solution for the quark propagator. 



ons are not necessarily massless in the chiral limit. For- 
mally, this problem is a higher-order effect in 1/N C but, 
as we will see below, it is quite severe. It can be avoided 
by performing a "strict l/A^-expansion", where only 
the Hartree gap equation is solved self-consistently and 
the NLO corrections are added perturbatively, discard- 
ing all higher-order terms. Then all diagrams contain 
only Hartree quark propagators and, hence, the inter- 
mediate RPA pions are massless in the chiral limit. As 
shown, e.g., in Refs. 

Hal EMI 

this pcrturbative treat- 
ment yields the correct results for the low-temperature 
behavior of the quark condensate and the pressure. How- 
ever, since the focus of the present paper is on the phase 
transition, which cannot be treated perturbatively, we 
stay with the self-consistent expansion scheme outlined 
above. We should then be alerted to the fact that the 
intermediate RPA mesons do not obey the chiral Ward 
identities. 



IV. NUMERICAL RESULTS 

In this section we present numerical solutions of the 
self-consistent gap equation at NLO and related quanti- 
ties. The main complication as compared to the leading- 
order problem arises from the fact that the NLO self- 
energy correction (last diagram in Fig. [5]) is non-local. 
Hence, unlike the Hartree self-energy, which only yields 
a constant shift in the mass, the self-energy is now energy 
and momentum dependent and consists of several terms 
with different Dirac structure. Assuming a homogeneous 
medium with even parity, the inverse propagator can be 
parameterized as 

S-\z,k) = j zC(z,\k\)-j-kA(z,\k\)-B(z,\k\), 

(20) 

where z is a complex energy variable. In vacuum, as a 
consequence of Lorentz invariance, the dressing functions 
A, B, and C are only functions of k 2 = z 1 — |fc| 2 , and the 
functions A and C are equal. In the medium, however, 
where we have a preferred frame, we have three indepen- 
dent functions, which depend on energy and momentum 
separately. 

A. Model parameters and computational details 

The integrals given in Sec. |H] are divergent and our 
model is only well-defined after specifying how to regu- 
larize them. Since the NJL model is non-renormalizablc, 
new cutoff parameters can appear at each loop order. For 
instance, even if we have regularized the the quark loops 
in the Hartree self-energy, Eq. ([3]), and the RPA polar- 
ization loop, Eq. (J8J), the loop over the meson momentum 
in I]' 1 ), Eq. (|T9"]) . is in general still divergent and needs 
to be regularized separately. 

In the following we regularize both, quark and meson 
loops, by sharp 3- momentum cutoffs A and Am, respec- 
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tively. This has the advantage that it keeps the numerical 
effort for the involved self-consistency problem as simple 
as possible. Moreover, it preserves the analytic structure 
in the complex energy plane. The obvious disadvantage is 
that the 3-momentum cutoffs violate the Lorentz covari- 
ance. For the moment, we take this as a minor problem, 
which could be improved on in future modifications of 
the model. 

The non-covariance of the regularization also makes 
it necessary to specify how external 3-momenta are dis- 
tributed to the propagators in a loop. In the RPA 
polarization functions, Eq. ([§])■ we distribute q equally 
to both quark propagators, whereas in the NLO meson 
loops, Eq. f(T9]) . we must attribute the entire external 3- 
momentum to the quark propagator in order to be con- 
sistent with the derivation of this diagram from the 
functional. 

In addition to the cutoffs, the model has two more 
parameters, namely the coupling constant and the bare 
quark mass. We take 

A = 664.3 MeV, GA 2 = 2.06, m = 5.0 MeV , 

which in leading order (Hartree / RPA) yield the empir- 
ical vacuum values for the pion mass m n — 135.0 MeV 
and pion decay constant f n = 92.4 as well as a quark con- 
densate of (to) 1 / 3 = -250.8 MeV [Hj]. This corresponds 
to a constituent mass van = 300 MeV. In principle, a re- 
fit of these parameters should be done at NLO. However, 
for the mostly explorative studies of the present paper, 
we keep them unchanged. The meson loop cutoff is set 
to A M = 500.0 MeV. 

The gap equation is solved iteratively starting with a 
Hartree-like ansatz for the quark propagator. The dress- 
ing functions are stored on a grid. In energy direction 
the grid is fixed through the Matsubara frequencies and 
in 3-momentum direction an equidistant grid space of 50 
MeV is chosen. Values in between the grid points are 
interpolated with cubic splines. The inverse propagators 
of the intermediate mesons are also stored on a grid. 



B. Dressing functions and intermediate meson 
propagators 

Results for the dressing functions A, B, and C arc 
displayed in Fig. [6] To a good approximation they can 
be taken to represent the dressing functions in vacuum, 
although for numerical reasons they have been calcu- 
lated at a temperature of 10 MeV and are therefore 
only given at discrete Matsubara frequencies. The results 
should be compared with the Hartree results, which are 
Bh = ran = 300 MeV for the present parameters and 
Ah = Ch = 1. As one can see, the NLO corrections lead 
to an overall reduction of the B function, whereas A and 
C are slightly enhanced. All dressing functions have in 
common that they are maximal at the lowest \oJ n \ and 
monotonously decrease with increasing \oj n \. (Note that 
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FIG. 6. Dressing functions C, A, and B of the quark propaga- 
tor at T — 10 MeV and /j, = as functions of the Matsubara 
frequency for different 3-momenta. 



the functions are symmetric in LJ n .) The same behavior 
can be observed with increasing 3-momentum. Asymp- 
totically, the NLO self-energy contribution vanishes, and 
the B function approaches a constant value coming from 
the Hartree diagram. For the same reason the C function 
converges to the trivial value of 1. 

In principle, we would expect the same behavior for the 
A function. In fact, because of Lorentz covariance, the 
functions A and C should be equal in vacuum. However, 
as a consequence of the non-covariant regularization, this 
turns out not to be the case. Moreover, when we consider 
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FIG. 7. Pion (upper panel) and sigma (lower panel) propa- 
gators at T = 10 MeV and fi — as functions of the Matsub- 
ara frequency for vanishing 3-momentum. The intermediate 
RPA propagators based on the self-consistent solutions of the 
NLO gap equation are compared with the propagators in the 
Hartree + RPA scheme. The numerical results at the dis- 
crete Matsubara frequencies are indicated by points whereas 
the lines correspond to fits according to Eq. (1221 1, 



the non-trivial parts A — 1 and C — 1, the symmetry 
violation is of the same order of magnitude as the physical 
effect. A closer inspection reveals that the A function is 
most strongly affected by the cutoff artifacts, because it 
is directly related to the 3-momentum (see Eq. (|2"0")l ). In 
particular, this explains the wrong asymptotic behavior 
of A. Hence, if we are interested in results which are 
sensitive to A — 1, an improved regularization scheme 
should be employed. For the quark condensate, which 
we discuss in Sec. IIVCI we expect that the situation is 
less problematic, as it is mainly influenced by the scalar 
function B (see Eqs. and (JM]) below). 

In Fig.[7]the intermediate RPA propagators in the pion 
(upper panel) and sigma (lower panel) channel channel 
at T = 10 MeV and \x — are displayed as functions 
of the Matsubara frequency for vanishing 3-momentum. 
The numerical results are indicated by points. For com- 
parison we also show the corresponding propagators in 
the Hartree + RPA scheme. As discussed at the end of 
Sec. IIIII the latter are constrained by chiral Ward iden- 
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FIG. 8. Intermediate sigma and pion propagators in the NLO 
scheme in the chiral limit at /i = and T — T c = 131 MeV 
as functions of the Matsubara frequency for vanishing 3- 
momentum. 



tities whereas the intermediate RPA propagators in the 
NLO scheme are not. In Fig. [7] this is reflected by the 
fact that the peak of the pion propagator at uj m = 
is strongly suppressed in the NLO scheme: In Hartree + 
RPA the peak is due to the relatively near-by pole at real 
energies, q° — — 135 MeV. The strong reduction of 
this peak in NLO is thus a hint for a considerably larger 
mass of the intermediate pion. 

We can estimate the meson masses by fitting the nu- 
merical points with a simple pole ansatz, 



Zm 



m M 



2G. 



(22) 



The constant — 2G has to be taken into account to get the 
correct asymptotic behavior (cf. Eq. (|4T|) below). These 
fits are indicated by the lines in Fig. [7] In Hartree + RPA 
we find 137 MeV, in good agreement with the true 

pole mass of 135 MeV. For the intermediate RPA pion in 
NLO, on the other hand, the fit yields m w ~ 340 MeV. 

In the sigma channel, the situation is less dramatic 
because the sigma meson is not a Goldstone boson. Here 
the pole fit yields m a « 670 MeV in Hartree + RPA and 
a slightly lower mass for the intermediate RPA sigma in 
NLO. However, these numbers should not be trusted too 
much, as they are based on rather far extrapolations from 
imaginary to real energies. (In fact, in the MEM analysis 
in the next section, we find masses which are 10 - 15% 
lower.) 

In Fig. |8] we show the behavior of the intermediate 
sigma and pion propagators in the chiral limit at the 
chiral restoration temperature, T = T c . As chiral sym- 
metry is restored, sigma and pion are now degenerate. 
However, even at T c the intermediate RPA mesons do 
not become massless. Fitting again the numerical points 
with Eq. (g2J, we find ra^ = m a w 270 MeV. This will 
be relevant for the discussion below. 
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FIG. 9. Temperature dependence of the chiral quark conden- 
sate in the chiral limit in Hartree approximation (dashed) and 
in NLO (solid) at fi = 0. 



C. The quark condensate 

The chiral quark condensate is a scalar quantity and 
can be calculated directly in Euclidean (Matsubara) 
space, 

r r1 3 k 

(qq)T = T Y,J ^s^(S(iu>„,k)), (23) 

where 



Tr(S(«w„, k)) 
-8N C 



B{iuj n , k) 



~uj 2 C 2 (iuj n , k) + k 2 A 2 (iuj n , k) + B 2 (iuj n , k) 



(24) 



cf. Eq. (|20jl. 

Our results for its temperature dependence in the chi- 
ral limit are displayed in Fig. [§1 both, in Hartree approx- 
imation (dashed) and in NLO (solid). We find that chiral 
symmetry is restored in a second-order phase transition 
in both cases. The critical temperature is decreased by 
the NLO corrections. 

A second-order phase transition is also what is ex- 
pected for two-flavor QCD in the chiral limit. At the 
critical temperature the fermionic degrees of freedom are 
suppressed due to their antiperiodicity. Therefore the 
phase transition is dominated by four bosonic degrees of 
freedom (three pions and the sigma meson) which are all 
massless at T c and one expects critical behavior according 
to the 0(4) universality class 18]. Since these arguments 
are only based on the symmetries and dimensionality of 
the system, the same should hold in any theory or model 
with the same conditions. Indeed, a second-oder phase 
transition with 0(4) critical exponents has been found in 
a renormalization group (RG) approach to the two-flavor 
quark meson model [23|, [24[ , and the N JL model should 
in principle exhibit a similar behavior. 



However, it is not a priori clear to what extent the 
universal behavior is spoiled by the truncation scheme. 
As we have seen, the massless bosonic degrees of freedom 
which are the basis of the universality arguments are not 
manifest in the gap equation, neither in Hartree approx- 
imation nor in NLO: In the Hartree self-energy, there is 
no back-reaction of the RPA mesons on the quark propa- 
gator, whereas the intermediate RPA mesons which enter 
the NLO gap equation are not massless. From this point 
of view, it seems not even guaranteed that the phase 
transition must be second order in these approximation 
schemes. To understand why this nevertheless should be 
expected, we can adopt the arguments of Ref . [25[ , where 
the chiral phase transition was investigated in a purely 
bosonic model at NLO and found to be second order as 
well: As we have discussed, at each order one can in 
principle construct mesonic correlators with the correct 
chiral behavior by applying Eq. (fT5]) to the ^-functional 
and iterating the resulting scattering kernel in the BSE. 
One could then employ these correlators to study critical 
exponents. Obviously, this can only work if the phase 
transition is second order. Thus, the gap equation must 
somehow "know" about the massless degrees of freedom 
even if they do not enter the equation explicitly. 

In this context the consistenc y o f the approximation 
scheme is crucial: In Refs. fl3l.1l5| a first-order phase 
transition was found in a simplified NLO scheme, which 
was suggested in Refs. 0, [Hj]. In that scheme, only lo- 
cal contributions to the quark self-energy are taken into 
account. Although one can formally show that there are 
massless Goldstone bosons in the chiral limit, the ap- 
proximation is not thermodynamically consistent as the 
gap equation cannot be derived from a thermodynamic 
potential. This suggests that a thermodynamically con- 
sistent treatment is important to find the correct order 
of the phase transition. 

Yet, even in a self-consistent and thermodynamically 
consistent truncation scheme, not all details of the critical 
behavior are necessarily reproduced correctly. It is well 
known that the Hartree gap equation, although correctly 
predicting a second-order phase transition, yields mean- 
field critical exponents. At NLO, there might be some 
improvement, but we should not expect to find exact 
0(4) behavior. Unfortunately, our numerical results are 
not precise enough to work this out quantitatively. 

The low-temperature behavior of the quark condensate 
is model independently given by chiral perturbation the- 
ory. For two quark flavors one finds in the chiral limit [26j 



rp2 rp4 

(qq)t = (qq)o ( 1 - ^77 



384/4 



(25) 



where (qq)o is the condensate at T = and / ff is the 
pion decay constant. This behavior is entirely due to the 
massless chiral pions. The T 2 term corresponds to their 
ideal gas contribution, while the T 4 term is due to p-wave 
7r — 7r interactions. 

Obviously, this behavior cannot be reproduced in the 
Hartree approximation, which includes no back-coupling 
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FIG. 10. Temperature dependence of the quark condensate 
for a bare quark mass mo = 5 MeV in Hartree approximation 
(dashed) and NLO (solid) at p = 0. 
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FIG. 11. Quark condensate in the chiral limit as a function 
of the chemical potential in Hartree approximation (dashed) 
and NLO (solid) at T = 0. 



of the mesons to the quark propagator. The change of 
the condensate is then exclusively triggered by thermal 
quarks and, hence, exponentially suppressed due to their 
mass. 

On the other hand, it was shown in Ref. [ll| that at 
least the T 2 term in Eq. (|25p is reproduced correctly, if 
the 1/N C corrections to the quark propagator are taken 
into account perturbatively. To be more precise, it was 
shown that the leading correction corresponds to that 
of an ideal gas of RPA pions. Since in the perturbativc 
approach the latter are built from Hartree quarks, they 
are massless in the chiral limit and, thus, lead to the 
correct low-temperature behavior. 

Accordingly, the T 2 -term in Eq. (|2"5|) cannot be re- 
produced in our fully self-consistent NLO scheme. As 
in the perturbative approach, there are corrections from 
the intermediate RPA mesons. However, because of their 
relatively large masses, their effect is exponentially sup- 
pressed. We are thus faced with the situation that the 
self-consistent scheme gives only a poor description of the 
low-temperature behavior but works well at the phase 
transition, while it is just the other way around in the 
perturbative approach. 

In Fig. [10] we show the temperature dependence of the 
quark condensate for the case of a non- vanishing bare 
quark mass. In this case, chiral symmetry is only ap- 
proximately restored in a crossover. Again, the NLO 
corrections lead to a reduction of the cross-over temper- 
ature relative to the Hartree result. 

The dependence of the quark condensate on the chem- 
ical potential at T = can be seen in Fig. Qj] for the 
chiral limit. The Hartree and NLO results look qualita- 
tively similar. In both cases chiral symmetry is restored 
in a first-order phase transition. Similar to the temper- 
ature behavior, the critical chemical potential is slightly 
lower in NLO. 



V. ANALYTIC CONTINUATION WITH THE 
MAXIMUM-ENTROPY-METHOD (MEM) 

Searching the analytic continuation of a propagator 
given at a discrete set of Matsubara frequencies is an 
ill-posed problem. In general, an infinite set of functions 
would provide an analytic continuation. Additional re- 
quirements on asymptotics and analytical structure pro- 
vide a unique solution [27| but this only helps if one can 
get analytic expressions for the functions. 

One way to attack this problem numerically is the 
Maximum-Entropy Method (MEM). The quantity to cal- 
culate is the spectral function p(uj) which is related to the 
propagator via the Lehmann representation (e.g. [H|), 



D(z) = 



dio p(ui) 
2n z — uj 



(26) 



and can be interpreted as a probability distribution. The 
method performs a x 2_n t to the discrete "data" D(iuj n ) 
and additionally requires minimum deviation from the 
so-called prior estimate, the "most probable" spectral 
function in absence of any data. The combination of both 
requirements finally leads to a unique solution. MEM 
has been successfully used in lattice QCD [29], [30| and 
Dyson- Schwinger calculations (3l| . 



A. The Method 

We largely follow the formalism described in Refs. [29], 
l3l| . The basis for MEM is B ayes' theorem for conditioned 
probability applied to the plausibility P[p\DH(m)] of a 
spectral function p under given data D and a prior knowl- 
edge H(m), 



P[p\DH(m)} 



P[D\pH{m)]P[p\H{7 
P[D\H{m)\ 



(27) 
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P[D\pH(m)] is the likelihood function, which indicates 
the plausibility of the data under the spectral function 
and the prior and P[p\H(m)] is the prior probability for 
the plausibility of p solely under the prior. P[D\H (m)] is 
a normalization factor independent of the spectral func- 
tion and can be dropped as the probabilities are normal- 
ized in the end. 

For the likelihood function a Gaussian distribution is 
assumed 



P[D\pH{m)\ 



1 



-L[P] 



with 



L[P] 



1 

~Nd 



E 



1A - D 
2\ai 



P I 2 



(28) 



(29) 



for Nd equidistant "data" with values Di and errors Oj 
and the corresponding values D? calculated from the 
given spectral function p using the Lchmann represen- 
tation, Eq. ([26)) . Zl is a normalization constant. Note 
that, by this assumption, P[D\pH (m)] does actually not 
depend on on the prior. 

The prior probability P[p\H(m)\ depends directly on a 
prior estimate m(w) for the spectral function. m(co) con- 
tains general information about spectral functions, es- 
pecially positivity, and is usually chosen as a constant 
function. With help of a scale factor a, P[p\H(m)] can 
be expressed as 



P[p\H(m)} 



daP[p\H(am)]P[a\H{m)} , (30) 



where 



P[p\H(am) 



Zs 



(31) 



with a normalization constant Zs and the Shannon- 
Jaynes entropy 



S[p] = dui f p(oj) — m(ijj) — p{uj) In 



m(oj) 



(32) 



This can be derived axiomatically by using general fea- 
tures of the entropy (locality, scale invariance etc.) or 
with help of the law of large numbers ("monkey argu- 
ment" ) [29( . Discretization and expanding the logarithm 
for small deviations of the spectral function from the 
prior yields 



S[p] « -2 Aw< iVPi ~ V^f ■ 

i 

Applying Eq. and Eq. $5T§ to Eq. (gTJ) gives 

1 



P[p\DH{am)] 



oQ[p\ 



(33) 



(34) 



with another normalization constant Z and the func- 
tional Q[p] — czS[p] — L[p]. Maximizing Q[p] gives the 
most probable spectral function, p(oj) = p a (oj), for given 
a. 

Finally the scale factor a has to be eliminated, which 
can be done in several ways. We use Bryan's method (3^ . 
which is applied in most cases. Here, the final spectral 
function Pmem{w) is obtained by averaging over a, 



Pmem(lo) = J Dpp{uj)P[p\DH{m)] 

OO 

w J dap a (u})P[a\DH(m)}, (35) 



with the probability factor 



P[a\DH{m)] cx J Dpe Q[p] 



exp 



2 V" V A fe 



+ Q[Po 



(36) 



and Xk being the eigenvalues of the matrix 



Mij = aAujiSij 



d 2 L 

'Pi imrrVPl 



1 dpidp 3 



(37) 



P=Pq 



In an intermediate step Laplace's rule (P[a\H(m)] = 
const.) [2§| has been applied. For calculating pmem{u) 
the probabilities Eq. ([36)) have to be normalized. 



B. Numerical implementation 

The iterative calculation of the quark propagator dis- 
cussed in Sec. HVI has been performed without error esti- 
mate, as the largest parts of the error are highly corre- 
lated and systematical. MEM requires data with noncor- 
related Gaussian errors (Eq. ([29])'). Therefore we assume 
a constant relative error of 10~ 4 all data points. This 
error underestimates the systematical errors but should 
be of the order the uncorrelated numerical errors. This 
choice seems useful as larger errors do not resolve the 
continuum contributions to the spectral functions while 
smaller errors lead to unphysical oscillations. 

We chose a constant prior of m(uj) = 10~ 3 McV -1 . 
This arbitrary choice can be justified as the result is quite 
insensitive, even when the prior is varied by six orders of 
magnitude, as illustrated in Fig. Q21 
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FIG. 12. MEM result for the NLO quark spectral function 
po{ui) in vacuum for vanishing 3-momentum, using different 
priors m(w) = const. 

C. Vacuum spectral functions 

Vacuurr0 results for the quark spectral function arc 
shown in Fig. [13] To be precise, since the spectral func- 
tion of spin-i fermions has a Dirac structure, we show 
the 0-component, 



1. 



pa = i Tr( 7 °p). 



(38) 

consider quarks with vanishing 3- 



Moreover, we 
momentum. 

The Hartree spectral function is displayed in the upper 
panel. Here we can directly compare the MEM result 
with the exact analytical solution, 



po(u) = n5(\u>\ - m H ) , 



(39) 



i.e., two delta peaks at w = ±300 MeV. Indeed, MEM 
yields two sharp peaks at the correct positions, although 
with a small width which is caused by the assumed nu- 
merical errors. 

The NLO spectral function (lower panel) has similar 
peaks, shifted to slightly lower energies, ui « ±270 MeV. 
In addition there are two broad bumps at higher values of 
\uj\. These structures can be related to the imaginary part 
of the NLO self-energy diagram (last diagram of Fig. [5]), 
i.e., to the continuum due to meson absorption or emis- 
sion processes on the quark. This continuum should be 
well separated from the quark mass peak and start at a 
threshold given by the sum of the quark mass and the in- 
termediate pion mass. According to our earlier estimate 
for m n , this would be at around 600 MeV and is more 
or less consistent with the MEM result. The details of 



2 As before, the "vacuum" results have been obtained at /i = and 
a temperature T = 10 MeV for numerical reasons. We expect 
the difference to real vacuum calculations at T = to be small. 
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FIG. 13. MEM results for the quark spectral function po(w) in 
Hartree approximation (upper panel) and NLO (lower panel) 
in vacuum for vanishing 3-momentum. 



the threshold region can, however, not be resolved. From 
Fig.[T2]we also see that the dip region, where we expect a 
vanishing spectral function, is most sensitive to the prior. 

Unlike in Hartree approximation, the MEM results 
for the spectral functions in the self-consistent NLO 
scheme cannot be confronted with direct calculations in 
Minkowski space, since the latter are highly non-trivial. 
However, such a comparison can be done for the quark 
propagator with perturbative NLO corrections to the self- 
energy. In that already mentioned, only the 
Hartree propagator Sh is calculated self-consistently and 
afterwards perturbative corrections are added to the self- 
energy, 



■^pert 



s 



H 



.(1) 
^pert i 



(40) 



where £p2-t i s strictly of the order 1/N C . In particular, it 
is entirely given in terms of Sh and meson propagators 
in the Hartree + RPA scheme. This makes its evaluation 
in Minkowski space possible. 

As discussed, e.g., in Ref. Spelt consists of two 
diagrams. The first one is local and, thus, only gives a 
constant real contribution to the mass. In the following 
comparison we will therefore neglect this term for sim- 
plicity. The second contribution is non-local and given 
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FIG. 14. Comparison between MEM result and direct 
Minkowski-space calculation of the quark spectral function 
with perturbative 1/N C corrections to the self-energy in vac- 
uum. 



by self-energy in the last diagram in Fig. [5] if all self- 
consistent quark propagators are replaced by Sh- In 
Minkowski space we make the additional approximation 
of neglecting the continuum parts and finite widths of the 
mass peaks of the RPA mesons. As the mass peaks are 
the dominant contributions to the spectrum and even for 
the unstable sigma meson the main peak is very sharp, 
this approximation should be quite accurate. 

A comparison with the MEM output of a correspond- 
ing perturbatively dressed propagator calculated in Eu- 
clidean space is shown in Fig. 1141 For numerical reasons 
the calculation in Minkowski space has been performed 
with quark propagators with a width of 1 MeV, which 
was added by hand. The main peak of the MEM re- 
sult and the Minkowski result are almost at the same 
energy. (The small difference could be due to the men- 
tioned approximations for the meson spectral function 
in Minkowski space.) The continuum contributions only 
show rough agreement in height and position but the 
shape differs considerably. The main reason for this is 
the dominance of the mass peak in the spectrum. The 
contina give only small contributions and it would re- 
quire data with much lower errors to become sensitive to 
the shape of the continuum. Furthermore the continuum 
structure is very complicated with sharp thresholds and 
peaks which are difficult to be reproduced by MEM. In 
fact, the high-energy thresholds and peaks are artifacts 
of the regularization and therefore not even physical. In 
principle, one could try to include these effects into the 
prior m{uS) as done in the lattice calculations of Ref. [33| . 

For applying MEM to meson propagators we use a sub- 
tracted dispersion relation, 



D M (z) = 



dui p(cj) 
2ir z — u 



-2G. 



(41) 
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FIG. 15. MEM results for the RPA pion (upper panel) and 
sigma (lower panel) spectral function in vacuum for vanishing 
3-momentum. The various curves correspond to the interme- 
diate RPA mesons in the NLO scheme for mo = 5 MeV and 
in the chiral limit, as well as in the Hartree + RPA scheme 
for m = 5 MeV. 



meson propagators. 

The spectral functions for the RPA meson propagators 
are displayed in Fig. [I3J The intermediate RPA mesons 
in the NLO scheme with finite current quark mass and 
in the chiral limit are compared with the mesons in the 
Hartree + RPA scheme. 

The latter can also be calculated directly in Minkowski 
space (cf. Sec. IV D[) . In the pion channel this yields 
a delta peak at u> = m v = 135 MeV and a contin- 
uum due to quark-antiquark decay above a threshold of 
2itih — 600 MeV. As one can see in the figure (upper 
panel, dotted line), these features are again qualitatively 
reproduced by MEM: We find a sharp peak at 135 MeV 
and a broad bump at higher energies, which can be iden- 
tified with the continuum. On the other hand, as before, 
the detailed threshold behavior cannot be reproduced. 

The spectral function of the intermediate RPA pion in 
NLO (solid line) looks qualitatively similar, but the mass 
peak is now located at u = 320 MeV. This is slightly 
lower than our estimate by the pole fit in Sec. IIVBI 
(340 MeV), but still very heavy. As we have said re- 
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FIG. 16. MEM result for the quark spectral functions p(ui) in 
Hartree approximation (upper panel) and NLO (lower panel) 
for vanishing 3-momentum at T — 200 MeV and \i = 0. The 
solid (dashed) lines correspond to an assumed relative numer- 
ical error of 10" 4 <->r,-6\ 



(10" b ). 



peatedly, this reflects the fact that the intermediate RPA 
mesons are not constrained by chiral Ward identities. In 
fact, even in the chiral limit (dashed line), we find a mass 
peak at 280 MeV in the MEM spectral function. 

In the sigma channel, the exact spectral function in 
Hartree + RPA has a mass peak slightly above the con- 
tinuum threshold. This is again qualitatively reproduced 
by MEM (lower panel, dotted line), where we find a max- 
imum at about 615 MeV. The NLO corrections to the 
sigma are rather small (solid line), the mass is lowered 
by about 20 MeV. In the chiral limit (dashed line), it is 
further reduced by 60 MeV, but the gross features of the 
spectral function remain unchanged. 



D. In- medium spectral functions 

Larger temperatures increase the inaccuracies signifi- 
cantly as the larger distance between the Matsubara fre- 
quencies (2ttT) provides only few data at low frequen- 
cies while high-frequency data mainly carry information 
about the asymptotics. 

As for larger temperature the quark mass decreases 



and therefore the main peaks in the spectral function 
come closer together, extrapolating the 0-component of 
the spectral function with MEM does not resolve the sin- 
gle peaks. Therefore it is reasonable to study other pro- 
jections of the spectral function to separate the peaks. 
For vanishing 3-momentum this can be achieved by the 
decomposition of the spectral function 



p{w) = p + (cj)L + 7o + p (lu)L 7q 



(42) 



with the projectors = i(l±7o). The components p^ 
correspond to the spectra of particle and antiparticles 
and are also positive definite, therefore MEM is applica- 
ble. In Hartree approximation the analytic result yields 



H J S(oj - m) 



p (w) =7T (1 



(43) 



(44) 



and the two main peaks are separated into the differ- 
ent channels. This decomposition has also been used in 
recent lattice studies at finite temperature [341 ] . 

MEM results for T = 200 MeV are shown in Fig. M 
While in Hartree approximation the quark and antiquark 
peaks are still well-separated the NLO peaks are much 
broader and overlap. This is an indicator for the thermal 
broadening of the quark spectral function in the medium. 

Applying MEM at finite chemical potential visualizes 
the shifted Fermi surface. The quark spectral function at 
p = 200 MeV is displayed in Fig. [T7J We obtain a quite 
asymmetric result: At positive energies we find a very 
sharp mass peak which is well separated from the con- 
tinuum at higher energies. At negative energies, on the 
other hand, both, the mass peak and the continuum are 
broader and the dip region in between is less pronounced. 

Although at finite chemical potential the spectral func- 
tion is indeed in general no longer symmetric, the main 
reason for the observed asymmetry is probably caused 
by the MEM procedure. This is due to the fact that the 
dispersion relation is now of the form 



S(iu n + p) 



dui 



p(w) 



2-7T iuj n + p — UJ 



(45) 



Therefore, the integral is most sensitive to details of the 
spectral function near the Fermi surface, where it gives 
the largest contribution. Hence, in the present example, 
the MEM procedure works best around oj = 200 MeV, 
which explains why the positive mass peak is much better 
resolved. 

In the mesonic sector, Hartree + RPA meson spectral 
functions can be calculated directly in Minkowski space 
without further assumptions. A comparison with MEM 
results in the pion channel at T = 100 MeV is displayed 
in Fig. [TH]for \q\ = (upper panel) and 100 MeV (lower 
panel). In the Minkowski-space calculations, a width of 
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FIG. 17. MEM result for the NLO quark spectral function 
po(oj) for vanishing 3- momentum at T — and fj, = 200 MeV. 



1 MeV was again added by hand to the propagators for 
numerical reasons. 

We find that the agreement of the MEM results with 
the direct calculations is rather poor. Only the posi- 
tion of the main peak fits while the shape of the spectral 
function looks totally different. The spacelike particle- 
hole branch in the spectrum which occurs at finite 3- 
momentum cannot be resolved in the MEM output. 

In order to shed some light on this problem, we show 
in Fig. Q1J] the Euclidean data which served as input for 
the MEM results in Fig. [THJ At low Matsubara frequen- 
cies the data points with different 3-momenta are slightly 
shifted against each other, which takes care of the differ- 
ent mass peak positions. But one cannot see a qualitative 
difference which could produce the spacelike continuum 
contributions which exist at \q\ — 100 MeV, but not at 
\q\ =0. Furthermore the complicated structure of the 
spectral functions is not visible in the Euclidean data 
and so this structure gets almost lost in the convolution 
with the Lehmann representation Eq. (|26|) . Instead the 
fit of a propagator with a single mass pole according to a 
delta peak in the spectrum already fits the data almost 
perfectly. This exemplifies how difficult it is to regain the 
spectral function from the Euclidean data. 



VI. SUMMARY 

We have studied the NJL model in next-to-leading or- 
der in a self-consistent l/A^-expansion. The chiral con- 
densate showes a second-order phase transition at finite 
temperature and vanishing chemical potential in agree- 
ment with expectations from the 0(4) universality class. 
In comparison to mean-field results the critical temper- 
ature is decreased. The low temperature behavior ex- 
pected from chiral perturbation theory cannot be repro- 
duced as the intermediate mesons which enter the NLO 
quark self-energy diagram are not restricted in terms of 
chiral symmetry and are very massive. At finite chemical 
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FIG. 18. Hartree + RPA pion spectral function in directly 
calculated in Minkowski space and extrapolated with MEM 
at T = 100 MeV for 3- momentum |g| = (upper panel) and 
\q\ = 100 MeV (lower panel). 



potential and vanishing temperature the phase transition 
is of first order like in mean-field approximation. The 
critical chemical potential is also slightly reduced. 

The Maximum-Entropy-Method has been used to cal- 
culate spectral functions for real energies from the Eu- 
clidean propagators. The mass peak of the particles can 
be reproduced quite well and for NLO quark propagators 
also a continuum contribution can be identified. The 
large mass of the intermediate pions is confirmed and 
approximately of order ~ 320 MeV while the masses of 
the quarks (~ 270 MeV) and intermediate sigma mesons 
(~ 600 MeV) are of the order of the mean-field values. 
Finite temperature increases the inaccuracies of MEM as 
fewer data at low energies are available. For perturba- 
tively dressed propagators the spectral function can be 
calculated directly in Minkowski space. A comparison 
with these results shows that the main peak is repro- 
duced well by MEM, but not the continua. 
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